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In two-dimensional (2D) inviscid incompressible flow, low background vorticity distribution 
accelerates intense vortices (clumps) to merge each other and to array in the symmetric 
pattern which is called "vortex crystals" ; they are observed in the experiments on pure 
electron plasma and the simulations of Euler fluid. Vortex merger is thought to be a result of 
negative "temperature" introduced by L. Onsager. Slight difference in the initial distribution 
from this leads to "vortex crystals" . We study these phenomena by examining iV-point 
vortex systems governed by the Hamilton equations of motion. First, we study a three-point 
vortex system without background distribution. It is known that a iV-point vortex system 
with boundary exhibits chaotic behavior for N > 3. In order to investigate the properties 
of the phase space structure of this three-point vortex system with circular boundary, we 
examine the Poincare plot of this system. Then we show that topology of the Poincare plot 
of this system drastically changes when the parameters, which are concerned with the sign 
of "temperature" , are varied. Next, we introduce a formula for energy spectrum of a Ap- 
point vortex system with circular boundary. Further, carrying out numerical computation, 
we reproduce a vortex crystal and a vortex merger in a few hundred point vortices system. 
We confirm that the energy of vortices is transferred from the clumps to the background in 
the course of vortex crystallization. In the vortex merging process, we numerically calculate 
the energy spectrum introduced above and confirm that it behaves as k~ a , (a as 2.2 — 2.8) 
at the region 10° < k < 10 1 after the merging. 

KEYWORDS: vortex, Poincare section, energy spectrum, vortex dynamics, vortex crystals, vor- 
tex merger, numerical simulation, negative temperature 
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1. Introduction 

It is well known that strongly magnetized pure electron plasma is approximately equivalent 
to 2D inviscid incompressible fluid. In a strong uniform magnetic field Bq = Bqz the drift 
velocity of the guiding center of plasma is given by v(x,y) = —V(j)(x,y) x z/Bq, where (j> is 
the electrostatic potential and z is a unit vector normal to the plane of the flow. Thus, this 
flow of the drift velocity is incompressible, i.e., V • v = 0. The vorticity of this velocity field 
is given by C( x >y) = V x v = zV 2 (J)/Bq. Using the Poisson equation, the vorticity of this 
velocity field is shown to be proportional to the electron density. Therefore, the vorticity is 
governed by the following equations, 



Here ((x,y) = (e/eoBo)n(x,y) ( n{x,y) is the electron density ) and \I/(x, y) = (eo/e)<fi(x,y) 
is the stream function, eq is the dielectric constant in vacuum. Thanks to the above relation, 
the behavior of vortices in incompressible Euler fluid is experimentally realized on non-neutral 
plasma. The vortex merging process and the motions of point vortices in vacuum or in the 
background vorticity distribution were experimentally observed. 1-3 

Theoretically, stability of point vortices array had been studied. 4,5 Furthermore equilib- 
rium states of iV-point vortex system had been also investigated. 6-9 The equilibrium states of 
this iV-point vortex system under various initial conditions were obtained by using statistical 
mechanical method. The theory using statistical mechanics predicted that the vortices of this 
system merge one after another to be a large vortex in equilibrium. This phenomenon occurs 
when the "temperature" of point vortices is negative. The positive "temperature" leads to 
merger between vortices of opposite sign, or with boundary. 8 

However, about ten years ago, a new phenomenon was observed in the experiment by 
Fine et al. w Under certain initial conditions, vortex merging is "cooled" and a few intense 
vortices (which is called "clump" in this paper to distinguish a point vortex) form a regular 
lattice in a low background vorticity distribution. This phenomenon which is called "vortex 
crystals" is the highlight of this field. This vortex crystals are investigated from various aspects. 
For example, a power law decay of the number of intense vortices of this system, namely 
N ~ = 0.2 — 0.7 was studied. 10 ' 11 The formation of the vortex crystals was predicted 

by two theories. One is a kinetic theory using the Langevin equation. 12 The other is a theory 
by the maximizing entropy principle. 13 

Moreover, this phenomenon, namely vortex crystal, which is experimentally observed 
was reproduced by numerical simulation. 14 This numerical result showed quantitatively good 
agreement with the experimental one. Time evolution of the number of intense vortices both 
in experimental flow and in numerical one has similar decay in time. 

In the recent experiment, 15 it is observed that three clumps arrayed on a straight line 
form an equilateral triangular lattice under certain conditions of the background. When the 
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background is very dense, merger immediately occurs instead of crystallization. In the very 
dilute background, chaotic motion of the three clumps which occurs without the background 
is not cooled by the background. Thus the crystallization does not occur and the three clumps 
in the dilute background do not exhibit a triangular configuration. Therefore the existence 
and the density of the background are supposed to be significant in the process of vortex 
crystals. 

In this paper, we report the results of our studies on these phenomena, namely the vor- 
tex crystals of three clumps, regarded as a iV-point vortex system. The vortex crystal was 
numerically reproduced by Schecter et al. lA However their simulation is based on the Euler 
equation (VIC simulation). 16 Therefore the details on the dynamics of vortices are not visible, 
namely how vortices in clumps and in background mix each other. Hence we treat the discrete 
iV-point vortex system to pursue the behavior of each point vortex. In the discrete iV-point 
vortex system, the detail of the dynamics such as mixing of vortices can be observed. In order 
to compare our simulation to the experiments on non-neutral plasma, the signs of circulation 
of vortices are all defined to be positive for our iV-point vortex system in this paper: in the 
experiments on non-neutral plasma, the sign of vorticities are all the same. 

The iV-point vortex system is described by the following Hamilton equations of motion 
for open boundary cases, 

H = ~4^X] ^2 r mTn^Og\z m - Z n \, (2) 

r dx m = dH_ r dym = dH 

m dt dy m ' m dt dx m 

Here z m = x m + iy m , i = (x m , y m ) and T m are the position and the circulation of the ru- 

th vortex, respectively. In this case, constants of motion are the total circulation ^2 m T m , the 
angular impulse I = J2 m r m{x m 2 + y m 2 ), the total energy E = -^EE( m ^) r ™ r n lo g|2m- 
z n \, and the center of vorticity Y, m T m z m / Yum r ™- 

For the system with circular boundary, the Hamiltonian is given by 

1 1 los R 

H = -^Y,Y1 r m r n io g |z m -z n | + — ^^r m r„iog|^ 2 -z m z„| — ^E r ™) 2 > ( 4 ) 

m n^m m n m 

where R is the radius of the circular boundary. The equations of motion are 

dx m \ Y n y m — y n \ Y n y m y n 

~df ~ ~ 2f 4^\z m -z n \ 2 ^4^\z m -z n '| 2 ' 

dj/m r ra X m X n \ T n X m X n . . 

dt ~ ^ 4vr \zm - z n \ 2 ^ 4tt \z m - z'\ 2 ' 



where 



/ _ R 2 i _ r2 i _ r2 / 7 n 
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The constants of motion for this system with circular boundary are the total circulation 
^2 m Fto) the angular impulse /, and the total energy E. 

The organization of this paper is as follows. In a circular domain, a system of more than 
or equal to three point vortices exhibits chaotic behavior. Therefore, as the most fundamental 
system of vortices, we investigate a three-point vortex system in §2 for the purpose of exam- 
ining how the system depend on the parameters (J, E), namely the angular impulse and the 
total energy. By plotting the Poincare map, we examine the dynamical behavior of this system 
in detail. We consider a parameter space which is spanned by two constants of motion (J, E). 
For various values of these two parameters, the Poincare plots are calculated. We confirm that 
the topology of the Poincare plot drastically changes when the values of the parameters are 
varied and that the property of the Poincare plot is influenced by the sign of "temperature" 
of point vortices introduced by L. Onsager. 8 In §3, we introduce an expression of the energy 
spectrum for this iV-point vortex system with circular boundary. This energy spectrum is 
used in order to understand the irreversible merging process in which a iV-point vortex sys- 
tem eventually grows into one large vortex, which is the result of negative "temperature" . In 
§4, we present the results of simulations. In a system of discrete point vortices, we reproduce 
vortex merger and a vortex crystal that is observed in the experiments on non-neutral plasma 
and in the simulations based on the 2D Euler equation. In §5, we summarize the results. 

2. Poincare Section 

In this section, we present the results of the study on a three-point vortex system. We 
investigate the topology of the Poincare plot for the three-point vortex system with 2D circular 
boundary under various initial conditions, which are determined by the parameters (I,E), 
namely the angular impulse and the total energy. It was confirmed experimentally and will be 
shown by our numerical simulation in the later section that the three clumps system exhibits 
a vortex crystal and vortex merger in the background vorticity distribution. Studying a three- 
point vortex system in vacuum is important in order to isolate the most fundamental properties 
from these phenomena, i.e., the vortex crystal and vortex merger: because it is supposed to 
be the backbone process of the three clumps system in the background vorticity distribution. 

It is known that in a system of point vortices without boundary chaotic behavior appears 
for more than or equal to four point vortices, and that with circular boundary it appears 
for more than or equal to three. Aref and Pomphrey 17 obtained the Poincare plot for a 
non-boundary system of four vortices and demonstrated that its dynamics is chaotic. Using 
their procedure, we calculate the Poincare plots for the dynamics of three point vortices with 
boundary. For the system with circular boundary, there are two constants of motion, namely 
the angular impulse / and the total energy E. We obtained the Poincare plots of this system 
for various values of these two constants. 

In the following, according to Aref and Pomphrey 17 we derive the Hamiltonian with two 
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degrees of freedom for our system by carrying out canonical transformations. For simplicity, 
we take the circulation T m = 1 for m = 1, 2 and 3. Then the Hamiltonian for this system is 
(The third term of eq.(4) is constant. So we ignore it in the following calculation.) 

^33 ^33 

H = ~ 5^ ^2 lo s( Zm ~ z n){z m ~ z n ) + — ^2 ^2 1 °S(R 2 ~ z m z n )(R 2 - z m z n ) (8) 



m=l n=l(^m) 



8vr 



m=l n=l 



Air 



\ogh(z\,Z2,zz), 



where 



h( Zl ,z 2 ,z 3 ) = (R 2 - \ Zl \ 2 ){R 2 - \z 2 \ 2 )(R 2 - \z 3 \ 2 ) 



\R 2 



- Z 1 Z 2 \ 2 \R 2 - Z 2 Z 3 \ 2 \R 2 - Z 3 Z X \ 2 

\zi - z 2 \ 2 \z 2 - z 3 \ 2 \z 3 - Zi\ 2 



• (9) 



The vortex coordinates {zi}i=i,2,3 are expanded in Fourier series 

75' 



(10) 



z 3 = ^(V2Jle w i + e" i4 f y^e^ 2 + e-^ V%h& 03 ). 

Here J„'s and ra 's are the action-angle variables. We further transform the variables in the 
following way. 

01 = |0?i -0 3 ), 

& = |(01-20 2 + 3 ), (11) 
03 = #2, 

V 

and 



a = Ji - j 3 , 
j 2 = Ji + j 3 , 
/ 3 = Ji + j 2 + j 3 - 

Substituting them into the Hamiltonian, we get the reduced Hamiltonian 

H I = H I (i 1 ,i 2 ,<h.,<h)- 

Here, since 2I 3 is just the angular impulse (i.e., the constant of motion, / = 2I 3 = ^ 



(12) 



m=l \'*m 

then 03 is a cyclic coordinate. Therefore it does not appear in the reduced Hamiltonian. 
Furthermore in order to get the Poincare section, we introduce canonical variables, 

R x = V(/i+/ 3 )/2/ 3 cos20i, 

< 

Pi = v / (/ 1 +/ 3 )/2/ 3 sin20 1 , 



(13) 

|2 ), 



(14) 



R2 = V(/ 3 -i2)/2J 3 sin202, 



P2 = V / U3-/2)/2/ 3 COs20 2 . 



(15) 
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The calculation of the Poincare plot is done by the following procedure. Since the reduced 
Hamiltonian eq.(13) is too complicated, time evolution is calculated by using the Hamiltonian 
eq.(8). We use the 4th order Runge-Kutta method for our numerical calculation. 10000 random 
initial points are given on the , Pi )-plane satisfying two constants of motion, namely the 
angular impulse and the total energy, and time evolution is done until these each orbit of 
initial points crosses the (Ri , Pi )-plane 50 times. The Poincare section is set on the plane 
with i?2 = 0, i?2 > 0. In this simulation, relative errors of both constants of motion, namely I 
and E, are 10~ 10 . 

First, we show the density of states as a function of the angular impulse I and the total 
energy E as Johnson did (Fig. 1). 18 A rigid line represents the values of the parameters at 
which three point vortices form an equilateral triangle configuration; the center of mass of the 
equilateral triangle is at the origin. The form of the rigid line is given by 



Here R is the radius of the circular boundary. The dotted lines indicate several contour levels 
of density of states. The innermost level is the highest (from outer to inner, 1.0 x 10~ 5 , 3.0 x 



10~ 5 , 1.0 x 10~ 4 , 3.0 x 10- 4 , 1.0 x 10~ 3 , 3.0 x 10~ 3 , 1.0 x 10~ 2 , 3.0 x 10~ 2 , 1.0 x 10 _1 (%) of all 



configurations, respectively). 

In the region for small / in Fig.l, the states with E < E T do not exist. The reason why 
the density in this region of parameters is zero is as follows. When vortices are close each 
other, the energy E increases. If vortices are close to boundary, the energy decreases. In small 
I range, because of the restriction of the angular impulse /, no vortices can be close to the 
boundary. Therefore, in the region of small I, the equilateral triangular configuration (rigid 
line) gives the minimum energy. This is actually observed in Fig.l. 

Figure 2 shows the Poincare plots for some parameters. The values of these parameters 
are indicated in Fig.3 and Table I. The density of the states reflects the area of the plots on 
the Poincare section. For example, we compare the areas of Fig.2(F) with that of Fig.2(G). 
Here we express the region of the plots in Fig. 2(F) as 2?( F ) on the complex plane Z = P\ + 



v/-LRi. Thus the region V (F) is expressed as P (F) = £>( F ) (0) U^(F) (1) U^(F) (2) U^(F) (3) , 
where £> (F ) (0) = {Z\0 < s < \Z\ < I}, and £> (F) (1) , £>( F ) {2) and £> (F) (3) represent three 



becomes thinner (Fig.2(G)). This tendency is observed for general values of I. 

When E = E r (hence the parameters (I, E) are on the rigid line in Fig.l), we can see from 
eqs. (10-15) that the plots on the Poincare section are lying on the circle of which radius is 



a/1/2 and at the origin. When the value of E changes from E r slightly for the same value of 

/, the region of the plots becomes concentric rings as Fig.2(A): if we write this two regions on 




(16) 



protuberances sprouted from V(p^ toward the origin. As the energy increases, three regions 
of the protuberances become thinner along the azimuthal direction. The region £> (F) (0) also 
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the complex plane Z = P\ + \f—\R\ as T>\ and T>2 , when E = E r these two region is expressed 
as T>\ = {Z\\Z\ = 0} and T>2 = {Z\\Z\ = a/1/2}. When E changes from E r , these regions 
change to T>\ = {Z\a\ < \Z\ < 02} and T>2 = {Z\b\ < \Z\ < 62}) where < a\ < a<i < b\ < 62- 
Therefore, the orbit in the phase space do not pass the origin on the (Pi, Pi)-plane. Actually, 
in Fig. 2, there are no points around the origin on the Poincare section. This means that the 
topology of the Poincare plot changes when the value of E crosses E r . As E increases, this 
two regions of the concentric rings are connected as seen in Fig. 2(B): these two regions T>\ 
and T>2 defined above change to T>\ = {Z\a\ < \Z\ < 02} and T>2 = {Z\b\ < \Z\ < 62}) 
where < a\ < b\ < CL2 < 62. As E increases more, the region T>i{jT>2 is scooped out: 
£>i U v 2 \ U™ 1 where £> (i) are the scooped region as seen in Fig. 2(C) (for m=6) and Fig. 
2(G) (for m=3). 

Furthermore, the property of Poincare plot is apparently changed when the sign of the gra- 
dient of the density of states dW/dE is changed to opposite sign. This gradient is proportional 
to the inverse temperature 1/T defined by "statistical mechanics". 

S = log W, 

1_ _ dS_ _ J_dW ( 17 ) 
T ~ dE ~ W~dE' 

where S is the entropy and E is the energy of the system. See Fig. 3. The values of the energy of 
points (A) and (B) are lower than the ridge of contour line (positive "temperature"), and (D) 
is around the ridge ("temperature" is infinite). The others are on the region upper than the 
ridge (negative "temperature"). When the "temperature" is negative, the region of Poincare 
plot is scooped out. This shows that the sign of "temperature" influences the behavior of the 
Poincare plot and the phase space structure. 

In the following, we take notice of the properties of the motion in the three-point vortex 
system with the parameter (I,E) = (100,2.25). Figure 4 shows the Poincare plot and the 
locus of three point vortices under different initial conditions on the real (x, y)-plane with 
(I,E) = (100,2.25). The plots of Figs.4(P-l) - 4(P-3) are generated from a single trajectory 
under different initial conditions. Figs. 4(0-1) - 4(0-3) show the motions of rotating three 
vortices. In the case of Fig.4(0-1), two vortices rotate each other around the origin apart 
from the other vortex which runs near the boundary. We will call the motion of this two 
vortices rotating each other in Fig.4(0-1) the "binary star motion". As seen in Fig.4(P-l), 
this motion of voritices plots tori on the (Pi, i?i)-plane. Tori on the Poincare section which 
does not correspond to the "binary star motion" of vortices are shown in Fig.4(P-2). The 
corresponding locus is depicted in Fig.4(0-2). Though the loci of these two type of motions in 
Figs. 4(0-1) and 4(0-2) are apparently different, these Poincare plots in Figs.4(P-l) and 4(P- 
2) are similar in its regular behavior on the Poincare section (Figs.4(P-l) and 4(P-2)). The 
other type of motions of vortices is shown in Figs.4(P-3) and 4(0-3). In this case, "chaotic sea" 
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emerges on the (Pi, i?i)-plane in Fig.4(P-3). The vortex orbits show cusps as in Fig. 4(0-3). As 
mentioned above, the motion of the three-point vortex system can be classified at least three 
types, namely the types of Figs. 4(0-1) - 4(0-3). The behavior of the Poincare plot for the 
motion of vortices enable us to distinguish the properties of vortices motion as in Figs.4(P-2) 
and 4(P-3): this difference is not distinguishable by only seeing the loci of three vortices in 
Figs.4(0-2) and 4(0-3). 

It is difficult to calculate the Poincare plot for the parameters in the region E < E r 
and in the region of large /. We suppose that there are two reasons why this calculation is 
difficult. The first reason is the complexity of the formula of the reduced Hamiltonian eq.(13). 
This complexity causes numerical errors. The second reason is the following. When E < E T , 
a vortex runs quite near the boundary (hence it runs very fast). In this situation we have 
to take the time step to be small. In actual numerical computation we fixed the time step. 
Therefore in that case numerical error becomes large in the calculation of the Poincare plot. 
This difficulty of numerical integration is a subject to be solved. A variable time step scheme 
would be needed. 

3. Energy Spectrum 

In a iV-point vortex system, when the temperature is negative, many like-signed vortices 
merge each other and tend to become a large vortex. 8 This indicates that the energy is 
transferred from small scale to large scale in two dimensional case. This phenomenon is known 
as the inverse cascade 19-21 in general two dimensional turbulence. 

Novikov gave a formula for the energy spectrum for a system of point vortices without 
boundary. 23 In the following, according to his procedure, we derive an explicit expression of 
the energy spectrum for a iV-point vortex system with circular boundary. 

The energy is defined by 



and £j(k) is the Fourier transform of the vorticity field. In this system, the velocity field is 
restricted to be within the circular domain. Thus we introduce an induced vortex sheet on 
the circular boundary whose radius is R. The vorticity field is given by 




(18) 



Here, v(k) is the Fourier transform of velocity field, 




(19) 




(20) 



rn 




(21) 



m 
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2irR\R-r m \ 2 

The effect of the induced vortex sheet is given by the second term in each line. Here r = 
(x,y), r m is the position of the m-th vortex, and R = Rr/r. *(r) is the stream function by 
vortices which is within the boundary, n is a unit vector perpendicular to the boundary, and 

• n = —d^/dr = vg. The derivation of the second term in eq.(22) is given in Appendix A. 

Thus the Fourier transform of the vorticity field uj(r) is given by, 

&(k) = — [ <jj(r) exp(ifc • r)dr 
2ir J 



5(r - r v 



R 2 ~ T r . 



2nR\R-r m \ 2 



5{r - R) 



This can be easily rewritten into 

^( fc ) = 2^£ r m I S(r - r m )exp(ife ■ r)dr - 

m 

1 ^ r R 2 — r 2 



R 2 -r 2 



2itR 



r=R 



2irR 

exp[ifc • (R - r m ) 
u r 12 

1 • ' m| 



exp(ifc • r)dr. 
exp(ifc • R 



(23) 
(24) 



r=R 



\R-r I 2 

I-"- 'm| 



dr 



■dr 



exp(ifc • r r 



^f m (fe)exp(ifc • r m ), 



where we define 



T m (^) — 27T^ r 



1 



2 „ 2 

' 7TJ 



cxp[ifc • {R-r m )\ 
\R-r m \ 2 



dr 



(25) 

(26) 
(27) 

(28) 



27T-R Jr=R 

In the limit A; — > 0, r m (fe) is vanished. Actually, the vorticity disappears in the outside of the 
circular domain. 



By definition, f m (fc) = f m (— k). Thus we have 

!~ (fe)| 2 = M*0P = J_^^f m (fe)f n (-fc)exp[ifc- (r m -r n )]. 

m n 

After some calculation, we have the formula for the energy spectrum, 

=^ I E r ™ 2 + 2 E E r m r n J (fc|r m - r n |) 

\ n m n<m 

m n 1=0 

oo j N 

££r m r n J>[j,(fci?)] 2 (^) cos^-^)]}, 



(29) 



(30) 



where <£> m is the angle between r m and the zero axis, and z m — r m e 1< ^ m . «//(x) is the /-th Bessel 
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function. e\ is denned as 



= < 



1 (1 = 0) 

2 (/>1). 



(31) 



Detail on the derivation is shown in Appendix B. Note that the first two terms in eq.(30) 
are just the formula derived by Novikov 23 and then the rest terms represent the effect of the 
circular boundary. 

4. Numerical Simulation 

In this section, we show the results of our numerical simulations on a vortex crystal and 
vortex merger. We confirm that the system of a few hundreds discrete point vortices is also 
able to exhibit a vortex crystal and a vortex merger. 

Non-neutral plasma experiments and 2D Euler equation numerical simulations found that 
by the existence of continuous low background vorticity distribution, in certain initial con- 
ditions intense vortices (clumps) are arrayed in a triangular lattice (vortex crystal). When 
the initial condition changes slightly to other distribution, clumps merge each other (vortex 
merger). This merger is the phenomenon seen when the "temperature" is negative. 

In our numerical calculations, the Hamilton equations of the iV-point vortex system are 
numerically integrated, not using the Euler equation. Numerical simulation is performed in 
the following procedure. We numerically integrate the equations of motion of point vortices 
with circular boundary (whose radius R=10), namely eqs.(5) and (6), by the 4th-order Runge- 
Kutta method. In our simulation, a relative error of the constants of motion is 10~ 7 — 10~ 8 for 
energy E and 10~ 10 — 10 -11 for angular impulse /. In the simulation using the Euler equation 
by Schecter et a/., 14 errors of both constants of motion are about 10~ 3 . 

4-1 Vortex crystals 

In the experiment by Sanpei et al, 15 it is found that three clumps are accelerated to form 
an equilateral triangular configuration by interaction with continuous low-level background 
vorticity distribution. Since three clumps make a unit cell of vortex crystals, it is regarded that 
this process is a fundamental process of forming vortex crystals. The existence of continuous 
low-level background vorticity distribution and the interaction between the background and 
the clumps are regarded as a key factor for making vortex crystals. 

We reproduce qualitatively a vortex crystal of three clumps in a system of point vortices. 
In Fig.5, our result is depicted. In our simulation, for the initial condition, each clump (densely 
distributed regions of point vortices) contains about 60 vortices, and the background distri- 
bution consist of about 180 vortices. The circulation of each vortex is 0.2. Three Clumps are 
arrayed on a straight-line with one clump at equal distances from others (Fig. 5 (a)) so as to be 
the same as the experiment by Sanpei et al. 15 The middle clump contains more point vortices 
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than others because of the overlap with background point vortices. 

Before seeing the result of vortex crystallization of our simulation, we define the property 
of clumps configuration, namely "arranging order". Labeling each clump by C{{i = 1,2 and 
3), rotating sequence of clumps in anti-clockwise is only (C1C2C3) or (C3C2C1). We call these 
rotating sequence the "arranging order". Vortex crystallization occurs as the following. At 
early steps, the clumps are at vertices of non-equilateral triangle rotating in the background 
point vortices (in Figs.5(a)-5(c)). Although the arranging order of the clumps do not change, 
the side length of the triangle of clumps changes. At later steps, the clumps form equilateral 
triangle and these rotate quasi-stationarily for a long time (we calculated until T = 200) (in 
Figs.5(d)-5(f)). 

In addition to the "arranging order", we introduce a signed symmetry parameter S. A 
symmetry parameter is defined as S = 12\^3A/l 2 by Sanpei et al. 15 A is the area and / is the 
peripheral length of the triangle whose vertices are the clumps. We define A as e z ■ (ri — r%) x 
(r2 — r^)/2. Here e z is a unit vector perpendicular to (x, y)-plane, and is the position of the 
i-th clump. An absolute value of S is maximized at 1 when the clumps form an equilateral 
triangle, and the sign of S is changed when the "arranging order" changes. In the following 
result of our simulation, the sign of S does not change because the "arranging order" does 
not change. 

We confirmed that the energy of clumps is transferred to the background vortices in the 
process of vortex crystallization (Fig. 6). The energy of three clumps decreases and the energy 
of the background increases while total energy of the system is conserved. In this process, the 
tendency of the energy of the background vortices is similar to that of the signed symmetry 
parameter (Fig. 7). Therefore, this means that if there is an outlet of energy of clumps, i.e., the 
background, three clumps tend to form a triangular lattice. However, it is necessary for the 
clumps not to occur merging. 

In this simulation, the arranging order of the clumps (namely, the sign of the signed sym- 
metry parameter S) does not change. However, in several parameters (the number of vortices, 
the circulation of background or clumps, and the ratio of these) this arranging order does 
change, i.e., the system shows the "binary star motion", and the system does not exhibit 
vortex crystals. In ref. 14, it is noticed that the ratio of the circulation of background to 
that of clumps influences the relaxation rate toward vortex crystals. According to the sim- 
ulation by Schecter et. al, 14 the relaxation rate toward vortex crystals is the fastest when 
Tbackground/rtotai is 0.2 - 0.4. In our case, r background /r t otai = 0.5. We suppose that the 
background distribution with the proper ratio to total circulation restrains the "binary star 
motion" to crystallize the clumps. 

It is observed that the point vortices of the clumps and the background mix each other 
(Fig.5). However, the number of vortices of each clump is conserved around 60. The vortices 
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of the clumps indicated by ▲ , T and ♦ are oozed out to the background. The vortices of the 
background indicated by + permeate to the clumps. 



4-2 Merging 

The background vortices causes not only crystallization but also merging. In this subsec- 
tion, we consider the iV-point vortex system as a model of vortex merging. Merging is highly 
important in a process of decay into vortex crystals, or into equilibrium states. Furthermore 
the power law behavior of the number of clumps (N ~ = 0.2 — 0.7) 10,11 in vortex crys- 

tals occur through the merging process of clumps. In merging process, we also calculate the 
energy spectrum of the system. In our simulation, each clump and the background contain 30 
point vortices. The circulation T m is set as T m = 1. For the initial configuration, the clumps 
are set at the vertices of equilateral triangle whose center of mass is at the origin. The energy 
spectrum is calculated by the formula derived in the previous section. 

We show the configuration of point vortices and the energy spectrum of this system in 
Fig. 8. The circular boundary is out of view (the area shown in Fig. 8 is the region [—6, 6] x 
[—6,6], the radius of boundary is R = 10). The behavior of peaks of the energy spectrum of 
the range 10° < k < 10 1 drastically changes in merging process. Before merging (in a lattice 
state, Figs. 8(a) and 8(b)), the peaks of the energy spectrum are oscillated with relatively large 
period in the range 10° < k < 10 1 . However, after merger (Figs. 8(c) and 8(d)) these become 
very sharp and wildly oscillates (i.e., plenty of peaks). 

The energy spectrum shows the power law E{k) ~ A; -2 ' 8 in the range 10° < k < 10 1 , and 
E{k) ~ /c~ L1 in the range k > 10 1 (Fig. 9). The origin of the power —1 region is the self energy 
(the first term of eq.(30)). The former power is near to the value of the prediction for the 
2D turbulence, 19-22 i.e., E(k) ~ k~ 3 . However, the other initial condition leads to different 
values, namely E(k) ~ k~ a ,a « 2.2 — 2.8. Moreover, the range of this power is changed for 
different distributions of vortices. The range 10° < k < 10 1 and the shape of spectrum in this 
range are an effect of the distribution of vortices in the merged clump. 

In order to confirm the effect of the background vorticity distribution, we calculate the 
system without the background of which configuration of vortices in the clumps are the same 
as the system with the background (Fig. 10). In this case, the clumps rotate preserving an 
equilateral triangular configuration for long time. Therefore this means that the background 
vortices merge clumps which are located separately so as not to merge without the background. 

5. Discussion and Conclusions 

In this paper, the vortex crystal and the vortex merger in the system of a few hundred point 
vortices are studied. The three-point vortex system without the background is investigated as 
a fundamental process of a three clumps system in the background vorticity distribution. 
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In §2, we studied the three-point vortex system with boundary. In order to investigate 
chaotic behavior of this system, the Poincare plot is calculated. The size of the plotted region 
on the section is closely related to the density of states: in the parameter where the density of 
state is low, the plotted area on the Poincare section is narrow and vice versa. The topology 
of the plotted area changes when the energy is varied across the value of the energy in which 
three point vortices are at the vertices of an equilateral triangular configuration. When point 
vortices are in an equilateral triangular configuration, the plots on the Poincare section are at 
the origin and on the circle of which radius is a/1/2. When the parameters change from this 
value, there are no plots at the origin on the section. We confirm that the sign of "temperature" 
influence the property of Poincare plot. When the "temperature" is negative, the region of 
Poincare plot is scooped out from inner of the region. The positive "temperature" leads to the 
Poincare plot in which there is no Tori. The property of the Poincare plot reflects the motion 
of point vortices. Tori or "chaotic sea" on this section have the relationship with the type of 
the motion of point vortices. The orbit that two point vortices rotate each other apart from 
the other point vortex draws tori on (Pi, i?i)-plane. The cusps are seen at the locus of the 
motion in which "chaotic sea" emerges on the section. The former indicates regular, and the 
latter shows chaotic behavior of point vortices. 

The energy spectrum of this system is introduced in §3. Using this result, we calculate the 
energy spectrum of the iV-point vortex system in merging process in §4. The energy spectrum 
in the range around 10° < k < 10 1 relaxes to the form of E{k) ~ k~ a , (a « 2.2 - 2.8). This 
range and the power are changed as the distribution of vortices is varied. In §4, the vortex 
crystal is reproduced in the system of point vortices. This shows that continuous distribution 
of vorticity is not necessary for vortex crystals. We find the energy transfer is occurred in 
the process of the vortex crystal of our simulation. The energy of clumps is disgorged to the 
background, and the clumps form triangle lattice. It is observed that the point vortices of 
the clumps and of the point vortices of the background are mixing each other. However, the 
number of voritices in each clump is conserved. This mixing effect of point vortices in the 
crystallization of clumps must be investigated. 

In our simulation, in the process of the vortex crystal, the arranging order of the clumps 
along the direction of rotation (namely, the sign of symmetry parameter S introduced at 
§4.1) does not change; however the arranging order of the point vortices are changed in a 
three-point vortex system as Fig.4(0-1). When this order of the clumps is permutated, vortex 
crystals does not occur in our simulations. This suggests a possibility that there are a relation 
between the crystallization of three clumps and the motion of three point vortices mentioned 
above. In a three-point vortex system, three vortices change the arranging order at some 
initial conditions and do not change at the other. We suppose that there are at least three 
factors for the crystallization of clumps, namely the parameters (/, E), the ratio of circulation 
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rbackground/Ttotai, and the ratio of the number of vortices ^background /^totai- Of course, the 
location of clumps also influences the crystallization. It must be investigated which parameter 
restrain the change of the order of rotating vortices. 
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Appendix A 



We want to calculate vg = —d^f/dr which is the second term of the vorticity field in 
eq.(21); the term of induced vortex sheet on the circular boundary. When the m-th vortex is 



at the position r T 



, 0), using the method of images, the stream function of this vortex 



at the position r = (r, 9) is 

Vm(r,0) = - — logr + — logr — log - — 

= - — log(r + r m - 2rr m cos 9) + — log(r + r m 
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+ 4^ l0g 



2\ 2 



r z + — -2r\ — cos 



i? 2 



r?rt . r 

2ir r m 



(A-l) 



where in the last line, the relation r m ' = R 2 /r rn is used, r' , r" and r m ' are, as given by Fig. 11, 
the distance between r and the m-th vortex, the distance between r and the image vortex of 
the m-th vortex, and the distance between the origin and the image vortex, respectively. 9 is 
the angle between r and r m . Therefore d^ m /dr is written as 



5*, 



2r — 2r m cos I 



2r - — cos I 



dr Air 
When r = R, eq.(A-2) becomes 
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2r m cos 9 



R 2 + r m 2 - 2Rr m cos 9 R 2 + r m 2 - 2Rr m cos ( 



2(R 2 



4ir R{R 2 + r m 2 - 2Rr m cos 9) 



= -r 



R — r^ 



m 2lrR\R-r m \ 2, 



(A-2) 

(A-3) 
(A-4) 
(A-5) 



where R = Rr/r. The sum about m of the above equation is the second term of eq.(22). 
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In this Appendix, the derivation of the energy spectrum E{k) is shown. The starting point 



is 



E(k) 



\v(k)\ 2 \u{k){ 



2 2k 2 
= 9I2 ^2^2^m(k)f n (-k)exp[lk ■ (r m - r„)], 



where, 



R 2 — r m 2 f exp(ifc • R) 



2ttR J r=R \R-r m \ 2 



dr ■ cxp(— ik ■ r m ) 



(B-l) 



(B-2) 



First, for the second term in eq.(B-2), we use the following relation 

exp(ifc • R) = Q ikRcos ( e -<t>m) _ gifc.Rsin(0-</> m -37r/2) 
00 

= Ji{kRy l{e ~^ m ~^' 2 \ 



(B-3) 



where 9 and (p m are the angles between r m and R, and between k and r m , respectively. 
Therefore we have 



r=R 



ex P (ifc-fi) dr = d9R Yr=-ooJi(kRy l{e -^ /2) 

\R-r m \ 2 Jo 



R 2 



1 + (it) 2 lt cos6 

'2vr 4l{0-4> m ) 



(B-4) 



1 ^ f^TT pU(V-Cp m ) 

- £ JKfcHXi) 1 / dfl -. (B-5) 

^z±^o ^ ! + -2^cos^ 



The imaginary part of the integral in the last line of eq.(B-5) vanishes because it is an odd 
function. The integral is evaluated as 



r 2w lie 
/ d9 1 

Jo i + ( r -f) 2 -r-?t 



Since r m /R < 1, the following relation can be used (\a\ < 1). 



cos ax ira 
■ ax 



1 — 2a cos x + a 2 



1 - a 



2 ' 



Thus, we have 

/ exp(ifc • R 

Jr=R \R~r, 



H 

CO 

J (kR) + 2 ^ MkR){i) 1 cos(/0 ri 



i=— 00 
2ttR 



R 2 ~ r 2 



1 m 

~R 



(B-7) 

(B-8) 
(B-9) 
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Finally eq.(B-2) is given by 

Z7T 



R 2 — r m 2 f exp(ifc • R) 



2ttR J r=R \R-r m \ 2 



dr ■ exp(— ifc • r m ) 



2vr m 



1 - ^2 eiJiikR) cos(l<p m ) ( ^ ) exp(-ifc • r m ) 



(B-10) 



where 



1 (Z = 0) 

2 (I > 1). 



Substituting this into |u)(fc)| 2 , we get 

; ( fc )! 2 = £X)f m (fc)f n (-fc)exp[i*- (r m -r n )] 



(27T) 

Here we define 

00 



expfifc • fr - r )} - V M - V W + V (m ' n) 
expire ^ m r n )\ 2—i 2—i 



(B-ll) 

(B-12) 
V . (B-13) 



} = 2^2 £ i' J i'( kRS ) cos ('Vn) ( — ^ ) exp(ifc-r m ), 



00 



Y^ m > = '^2 £l j l (kR)cos{l(f) m ) ( ^ ) exp(-ifc-r n ), 



1,1' l=0l'=0 

Integrating over angles ip of k, the energy spectrum can be written as 
r2 * \v(k)\ 2 1 ^ 



1=0 

00 00 



R 



(B-14) 



R \ R 
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1 /* 
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where, 
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(B-16) 



( A ) = 2 ( 27r )2 fc / d^^^r m r„exp[ife- (r m -r n )}, 

V / JO m n 

{B) = ~2{2^yki d ¥'EE r "» r »E e ' J '( fcl2 ) cos ( Z ^)(^J exp(-ife-r n ), (B-17) 
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1 [ 27T °° ( -\r \ V 

(C) = - 2(27T) 2 k / dip^2^2r m r n ^2ei,Ji,(kR)cos(l'(f) n ) ( — - j exp(ifc • r m ), 

(B-18) 
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2(2vr) 2 fc j 



00 00 
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x <^ 'Y^^ei£ V Ji(kR) cos(l4> m )Ji'(kR) cos(l'(f> n ) (^tt) l-jr 

{ 1=0 l'=0 \ / \ 

We calculate each term as follows. 

In order to calculate the first term (A), suppose that the angle between k and (r m — r n ) 

is ip. We can expand exp[ifc • (r m — r n )\ in the series of the Bessel functions, 

00 

exp[ifc • (r m - r n )\ = exp[ifc|r m - r n \ cosip] = ^ J^(fc|r m - r^De 1 ^ 7 ^. 

i=— 00 
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where ip = x — V and X is the angle of r m — r n . Then, the first term (A) is written 
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The second term (B) of eq.(B-15) is 
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(5) 
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Suppose that ip and v? n are the angle components of k and r n respectively. We have 

00 
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l=— 00 
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cos 



s (/0„) = I ^iv-vn) + e -u(v-iPn) J j (b-25) 
where n = - (/? n . Substituting these into the sum over I in eq.(B-23), we obtain 
J2eiJi(kR) cos(l(f> m ) ( ^ J exp(-ifc • r n ) 
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1 OO , . \ I oo 

= ^Y, ei yjr) Jl ( kR } E ^'(^ n )e- u '(^-^ + i){e"^-^)+e- y ^-^)}. (B-26) 

1=0 ^ ' l'=-oo 

Integrating over the angle ip, the sum over I' in eq.(B-26) is written as 

l>=-oo J ° 1 J 

= Ji{kr n ){-i) [ • 2cos[/(^ m - tp n )] • 2tt, (B-27) 

where the integral over <p of e ll P^ l ~ 1 ') is zero when I / I' and is 2ir when I = I'. Then we can 
write the second term (B) as 

^""^rifcSE^-E 6 '© J,(fei2) J,(fcr n ) ■ costi^m - Vn)]- (B-28) 

m n 1=0 

For the third term (C), a similar calculation to the second term (B) can be applied. Then 
we have 

1 00 1 
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The last term (Z>) is 
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Since I and Z' are both integers ( I, I' > ), then we have 
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(otherwise). 
Thus the fourth term (D) is given by 
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Substituting these four terms into eq.(B-15), the formula for the energy spectrum is given by 

E( k ) =4^ { E F " 2 + 2 E E r™r n J (A;|r m - r n \) 

\ n m n<m 

^EE^-E^l) Ji(kR)Ji(kr n ) ■ coa[l{<p m - <p n )] (B-33) 

m n /=0 

+EE r ™ r «E £ ' i j k^)] 2 (^) 
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Fig. 1. The contour plot of the density of states as a function of / and E calculated from 5.0 x 10 7 
random configurations of three vortices. A rigid line in this picture shows the values of parameters 
at equilateral triangle configuration of which center of mass is at the origin. The dotted lines 
indicate the contour levels (these are from outer to inner, 1.0 x 10~ 5 , 3.0 x 10~ 5 , 1.0 x 10~ 4 , 3.0 x 
10~ 4 , 1.0 xlO" 3 , 3.0 x 10~ 3 , 1.0 x 10~ 2 , 3.0 x 10~ 2 , l.Ox 10 _1 (%) configurations of all, respectivelly.). 
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Table I. The values of the parameters (I,E) for the Poincare plots shown in Fig. 2. E T indicates the 
value of the energy for the equilateral triangle configuration. 
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Fig. 2. The Poincare plots on the (Pi, i?i)-plane for various values of parameters. The values of 
parameters are summarized in Table I. (E-2) is the enlargement of (E). 

21/30 



J. Phys. Soc. Jpn. 



Full Paper 









1 


2.8 


I (Q 
v * 






2.6 






(G) 




\<B) 


•(E) 


• 


2.4 








2.2 




(D) 


(F) 






1 





50 100 150 200 

I 

Fig. 3. The function E?(I) in which three vortices form an equilateral triangle configuration (rigid 
line) and the values of the parameters for Poincare plots in Figs.2(A)-2(G). (A)-(G) correspond 
to the sections (A)-(G) in Fig. 2, respectively. 
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Fig. 4. The Poincare plot of a single orbit in phase space (P-1,2,3) and the locus of vortices on real 
plane (0-1,2,3) under different initial condition with (I, E) = (100,2.25). (P-1,2,3) are the parts 
of Fig.2(D). 
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Fig. 5. The vortex crystal of point vortices. Each clump contains 60 vortices. The background consists 
of 180 vortices. The circulation of each point vortex is 0.2. 



24/30 



J. Phys. Soc. Jpn. 



Full Paper 




Fig. 6. The energy of three clumps, the energy of the background vortices and the total energy of 
the system. The energy of the clumps tends to transfer toward the background. In later steps (the 
time crystallization occur), the energy transfer saturates. 
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Fig. 7. The signed symmetry parameter and the energy of the background vorticity. These show 
similar tendency that when the clumps form triangular lattice these value are growing. 
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Fig. 8. The configurations of vortices and the energy spectrum in the merging process. Each clump 
and the background consist of 30 point vortices. Around T = 2.0, a vortex merger occurs. 
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Fig. 10. The clumps without background vorticity distribution. Each clump contains 30 vortices. The 
initial configuration of vortices of clumps are same as Fig. 8. 
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Fig. 11. The image vortex 
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